Orientational Defects in Ice Ih: An Interpretation of Electrical Conductivity 
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We present a first-principles study of the structure and energetics of Bjerrum defects in ice Ih 
and compare the results to experimental electrical conductivity data. While the DFT result for 
the activation energy is in good agreement with experiment, we find that its two components have 
quite different values. Aside from providing new insight into the fundamental parameters of the 
microscopic electrical theory of ice, our results suggest the activity of traps in doped ice in the 
temperature regime typically assumed to be controlled by the free migration of L defects. 

PACS numbers: 61.72.Bb, 72.80.-r, 31.15.Ar 



While the isolated water molecule is one of the simplest 
in Nature, the condensed phases of H2O reveal many 
complex features that still elude complete understand- 
ing Q, . An example is one of the most abundant crys- 
talline solids on Earth, the proton-disordered hexagonal 
ice Ih, for which several aspects of structure-properties 
relationship have yet to be clarified Q. 

An important issue concerns the role of crystal defects 
in the peculiar electrical properties of ice Ih. When an 
electric field is applied to an ice specimen, it becomes 
polarized by the thermally activated reorientation of the 
molecular dipoles. To explain the molecular origin of 
this phenomenon, Bjerrum 0] postulated the existence 
of orientational defects that represent local disruptions 
of the hydrogen-bond network of ice Ih. While the con- 
ceptual picture of these Bjerrum defects is now well es- 
tablished, [3| a quantitative understanding of their struc- 
ture and energetics is still lacking, rendering a direct in- 
terpretation of experimental electrical conductivity data 
difficult. 

Recent atomistic studies ~>. (i. ?} based on empirical po- 
tentials have provided important qualitative insight into 
the structure and dynamics of Bjerrum defects, but have 
not yet attempted to make direct contact with experi- 
mental conductivity data. Furthermore, the few ab initio 
studies 0, in the literature involved clusters that are 
too small to reliably capture the properties of a defect 
embedded in bulk crystal. In this Letter, we present an 
ab initio study of the structure and energetics of Bjerrum 
defects in ice Ih using a large supercell. Moreover, based 
on Jaccard's defect-based microscopic electrical theory of 
ice 0, Ho) , we interpret the results in terms of experimen- 
tal electrical conductivity data for doped ice Ih. 

Fig. [IJi) provides a schematic picture of the forma- 
tion of Bjerrum defects in ice Ih. Defect-free ice Ih, 
in which each molecule is hydrogen-bonded to its four 
tetrahedrally positioned neighbors, obeys Pauling's two 
ice rules (i) each molecule offers/accepts two hydro- 



gens to/from two neighboring molecules, and (ii) there 
is precisely one hydrogen between each nearest-neighbor 
pair of oxygens. The proton-disordered character implies 
that there is no long-range order in the orientations of the 
H2O molecules or hydrogen bonds. The Bjerrum defect 
pair is a violation of ice rule (ii): it is obtained by the 
indicated rotation of molecule 1, after which there are 
two hydrogens between molecules 1 and 2 (D defect) and 
none between molecules 1 and 3 (L defect). After the 
initial formation, the pair can be further separated by 
successive molecular rotations, eventually creating a pair 
of independent, or free, D and L defects. Their motion 
through the crystal provokes the rotation of the water 
molecules along their paths, providing a mechanism for 
electrical polarization. 

In both pure as well as doped ice, the conductivity 
(Tdl due to Bjerrum defects is essentially controlled by 
L defects According to Jaccard's electrical theory of 
ice 01 it takes the form ctdl = <Zdl«lMl, where (Jdl 
is the effective charge carried by the Bjerrum defects, 
«l is the concentration of free L defects and is their 
mobility. Its temperature dependence is described by Q 



ctdl(T) 



T 



exp(-E/k B T), 



(1) 



where T is the absolute temperature, k B is Boltzmann's 
constant and E is a characteristic activation energy. In 
pure ice, the latter is given by Q 



E — tj-Edl 



E 



Lm j 



(2) 



where -Edl is the formation energy of an independent 
defect pair and Ei, m is the migration energy of a free L 
defect. 

While the activation energy can be determined experi- 
mentally by measuring ctdl as a function of T and adjust- 
ing it to Eq. JJJ, it is not possible to directly isolate the 
components Edl and E^ m . To achieve this, additional 
conductivity measurements need to be carried out on 



2 



doped ice samples, |3j in which an extrinsic concentration 
of L defects is injected. This indirect procedure, how- 
ever, is subject to large uncertainties and has revealed 
incompatibilities between different experiments 

Here we use a density-functional theory (DFT) 0] 
approach to explicitly compute the formation and mi- 
gration energies -Edl and E^ m . Our calculations are 
performed using the 96-molecule supercell labeled 3 x 
2 x 2 in Ref. 12. The calculations are executed us- 
ing the VASP package 0, Q] using the Perdew- 
Wang 91 generalized-gradient approximation ^l| and the 
projector-augmented-wave |l5| approach. Brillouin-zone 
sampling was limited to the T-point and we use a plane 
wave cut-off of E cut = 700 eV. The effects of spurious 
image dipole interactions were evaluated ^| and found 
to be negligible for all investigated structures. 

First, we relax the defect-free crystal supercell, allow- 
ing both the atomic and supercell coordinates to relax 
at zero stress. The resulting hexagonal lattice parame- 
ters a = 4.383 A and c = 7.16 A, are ~ 2 % below the 
experimental values measured at T — 10 K [3j. The aver- 
age intramolecular oxygen- hydrogen separation of 1.01 A 
is in excellent agreement with the experimental value of 
- 1.006 - 1.0081. As appears typical of DFT calcula- 
tions on ice Ih 17] our calculations slightly overestimate 
the sublimation energy of 0.69 eV compared to the ex- 
perimental value 0.61 eV Q. 

Next, we create an "embryonic" Bjerrum defect pair 
according to Fig. and relax it at constant volume. 
The resulting structure is shown in Fig. ^p) and is qual- 
itatively similar to those observed in a recent molecular 
dynamics (MD) study 0- Compared to the initial ge- 
ometry in Fig.^>), one can no longer recognize a D de- 
fect in the sense of its description in Fig. QJi) due to the 
large electrostatic repulsion between the two hydrogen 
atoms 0, 0, ■ The total-energy part of the formation 
energy of this structure is found to be 0.55 eV. To ac- 
count for zero-point contributions, known to be relevant 
in ice we also evaluate the change in the local in- 
ter and intramolecular vibrational modes with respect to 
the defect-free crystal for the molecules in the vicinity of 
the defect, using the local harmonic approximation [18j. 
Overall, the zero-point contribution lowers the formation 
energy by about 10% to 0.50 eV. 

To estimate the formation energy -Edl of an inde- 
pendent pair of Bjerrum defects, we move the L defect 
through the crystal by a series of molecular rotations, fol- 
lowed by structural relaxation. Assuming that the zero- 
point contribution is the same as for the embryonic defect 
pair, the formation energy as a function of the number of 
molecular rotations is shown in Fig. [5] It quickly reaches 
a plateau value around 1.1 eV, subject to small fluctu- 
ations due to the disordered character of the hydrogen- 
bond network. However, a more accurate estimate of 
Eul requires an analysis of its asymptotic behavior as 
a function of distance between the defect pair near the 



plateau value. Because of their effective charges ±<?dl, 
this behavior is expected to be of the form 

E iorm (r) = E DL - - g ° L , (3) 

4 7T60 Coo?" 

where r is the distance between the point charges ±<?dl, 
and Coo is the high-frequency dielectric constant of ice 
Ih For each relaxed Bjerrum defect pair we po- 

sition a point-charge +quh on the dangling proton of 
the D defect and center a charge — quh midway between 
the two oxygen atoms of the L defect, after which we 
adjust the formation energies of the defect pairs sep- 
arated by more than 3 molecular rotations to Eq. 10 • 
The results, shown in the inset of Fig. |5J are consistent 
with 1 /r behavior and the intercept gives an asymptotic 
value Eul = (1.153 ± 0.04) eV. In addition, the slope 
of the fit provides an estimate for the effective charge 
which, using the experimental value = 3.2 3], gives 
9dl = (0.34 ± 0.07) e, in good agreement with the exper- 
imental value (/dl = 0.38 e Q. 

Having determined £"dl we now compute the migra- 
tion energy barrier E^ m . Given the disorder in the 
hydrogen-bond network, this barrier is expected to fluc- 
tuate depending on the local environment of the L defect. 
For this purpose, we computed 6 distinct barriers at 3 dif- 
ferent L-defect sites (13, 14 and 17 molecular rotations, 
cf. Fig. [21 by starting from the relaxed D-L configura- 
tion and rotating one of the 2 molecules hosting the L 
defect as shown in Fig. Et): constraining only the rotat- 
ing proton while allowing full relaxation of all other coor- 
dinates. This yields transition states of the type shown 
in Fig. Ob), in which the angle ZABC is essentially bi- 
sected by the coplanar rotating OH bond. The resulting 
energy barriers (neglecting zero-point effects) vary in the 
range E^ m ~ 0.10 — 0.14 eV, indicating the significance 
of the disorder but an overall high mobility of free L de- 
fects. These results are consistent with the recent MD 
findings obtained at T = 230 K for which, assuming 
an attempt frequency vq = 24 THz typical for librational 
modes in ice Ih the present barriers give an average 
migration time ranging between 4 and 35 ps. 

The experimental estimates for the formation and mi- 
gration energies are obtained from conductivity mea- 
surements on ice doped with a substitutional concentra- 
tion of HF US IH or HC1 [H molecules. Since each 
molecule has only one proton, they introduce an extrin- 
sic, temperature-independent concentration of L defects. 
Measurements of the conductivity odl as a function of 
T then typically yield an Arrhenius-type plot with 3 dis- 
tinct activation energies, as shown in Fig.^Q- 

The high-temperature regime / is believed to be con- 
trolled by intrinsic behavior, characterized by the activa- 
tion energy Eq. J5J). The lower-temperature regimes II 
and III are assumed to be dominated by the extrinsic L 
defects. In the former, the temperature dependence of 
the conductivity has been attributed to the free motion 
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of extrinsic L defects so that £"11 = E^ m . In regime III 
the temperature is so low that the extrinsic L defects are 
not completely dissociated from their dopant molecules 
so that the activation energy involves an additional dis- 
sociation energy, giving E m = ^E diss + E Lm [jj. 

Considering the experimental results 0, 0, re- 
produced in Table 01 the activation energy for regime I is 
quite well established, even for different dopant species, 
showing a dispersion of less than 0.05 eV among the dif- 
ferent experiments. The values for the individual compo- 
nents -Edl and i?Lm, however, show considerably larger 
deviations. The migration barrier values vary between 
0.19 eV and 0.315 eV, leading to formation energies -E-dl 
ranging from 0.66 eV to 0.79 eV. 

Comparing these to our DFT calculations (cf. Table 
[IJ, we notice that both components of the activation en- 
ergy deviate significantly from the experimental values. 
The DFT result for Edl is more than 46 % larger than 
the largest experimental value, whereas E\^ m is about 
37 % lower than the lowest experimental estimate. In 
this light, it is quite striking that, despite the discrep- 
ancies for the individual components, the DFT estimate 
for the net activation energy in Eq. J5J agrees quite well 
with the experimental results, deviating about 10 % from 
the highest experimental estimate. This seems to demon- 
strate the difficulty involved in the interpretation of con- 
ductivity experiments in doped ice samples under condi- 
tions not controlled by intrinsic properties. Specifically, 
the fact that the DFT estimate for the migration barrier 
is systematically and significantly lower than all experi- 
mental estimates indicates that, as suggested in Ref. 0, 
the regime interpreted as being controlled by free extrin- 
sic L defects, may in fact involve the activity of traps 
that obstruct their motion, leading to the higher effec- 
tive migration barriers deduced experimentally 

In summary, we have conducted a first-principles study 
of the energetics of Bjerrum defects in ice Ih and com- 
pared the results to experimental electrical conductivity 
for doped ice samples. The results provide new insight 
into the parameters in Jaccards microscopic electrical 
theory of ice. While the DFT value for the net acti- 
vation energy is in good agreement with experiment, we 
find that its two components have quite different values 
from those inferred from experiment. In particular, our 
results predict a migration barrier for L-defect motion 
that is significantly lower than the lowest experimental 
estimate, hinting at the presence of traps in the regime 
usually interpreted as being controlled by the free migra- 
tion of extrinsic L defects. 
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TABLE I: Activation energies and their components Eul and 
i?Lm determined from doped ice experiments and the present 
DFT calculations. 





E 




-Bdl 


Ref. 10 (HF) 


0.575 


0.235 


0.68 


Ref. 19 (HF) 


0.625 


0.315 


0.664 


Ref. 20 (HC1) 


0.585 


0.190 


0.79 


This work 


0.696 ±0.04 


0.120±0.02 


1.153±0.04 
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FIG. 1: a) Schematic picture of the formation of a Bjerrum 
defect pair from the defect-free ice Ih structure. The defect 
pair is created by the indicated rotation of molecule 1. b) 
Relaxed DFT structure obtained from the configuration in 
panel a). 
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FIG. 2: Formation energy of the Bjerrum defect pair as 
a function of the number of molecular rotations separating 
them. Dashed line is to guide the eye. Inset shows least- 
squares fit of the formation energies of the defect pairs sepa- 
rated by more than 3 molecular rotations to Eq. • 
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FIG. 3: Identification of transition state for the migration 
of an independent L defect, a) Equilibrium structure of the 
independent L defect, b) Typical transition state. 
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FIG. 4: Schematic representation of a characteristic experi- 
mental Arrhenius plot of the conductivity <tdl as a function 
of temperature in a doped ice Ih sample 3] . 



